home *** CD-ROM | disk | FTP | other *** search
/ Tech Arsenal 1 / Tech Arsenal (Arsenal Computer).ISO / tek-20 / usat92.zip / USATQB.BAS < prev   
BASIC Source File  |  1991-12-07  |  26KB  |  885 lines
  1. ' AMSAT ORBITAL PREDICTION PROGRAM de W3IWI - May 1980
  2. ' COPYRIGHT 1980 by Dr. Thomas A. Clark, W3IWI
  3. '                   6388 Guilford Road
  4. '                   Clarksville, MD. 21029
  5.  
  6. ' REVISED & MODIFIED FOR IBM-PC by R. D. Welch, W0SL - May 1983
  7. '                                  908 Dutch Mill Dr.
  8. '                                  Manchester, Mo. 63011
  9.  
  10. ' ENHANCED AND DEBUGGED BY Ing. H.F.STRASSER, OE1HSI - JAN. 1985
  11. '                          A 1238
  12. '                          VIENNA, AUSTRIA
  13.  
  14. ' Modified to use the NASA/NORAD two-line element set format (as distributed
  15. ' by Dr. Kelso of the Air Force Institute of Technology) and revised for
  16. ' QBASIC by Antonio Querubin, AH6BW - Dec. 1991
  17. '    Internet:  tony@mpg.phys.hawaii.edu
  18. '    AMPRnet:  ah6bw@uhm.ampr.org
  19. '    PBBS:  ah6bw@ah6bw.hi.usa.oc
  20. '    BITnet:  querubin@uhunix
  21.  
  22. ' Permission granted for non-commercial use providing
  23. ' credit is given to the author(s), AMSAT and ORBIT Magazine.
  24.  
  25. DEFDBL A-Z
  26.  
  27. DECLARE FUNCTION gets$ ()
  28. DECLARE SUB getkey (k$)
  29. DECLARE SUB timezone (d0$, t0$, d1$, t1$, tzcor!)
  30. DECLARE FUNCTION juliantodate$ (day%, year%)
  31. DECLARE FUNCTION datetojulian% (d$)
  32. DECLARE FUNCTION trim$ (i$)
  33. DECLARE FUNCTION readtle% (tlefile$, name$(), epochyear%(), epochjulianday#(), inclination#(), TLERAAN#(), eccentricity#(), TLEArgOfPerigee#(), TLEMA#(), MeanMotion#(), TLERevolution&())
  34. DECLARE SUB readgrnd ()
  35. DECLARE SUB plotloc (longitude!, latitude!, xcoord%, ycoord%)
  36. DECLARE FUNCTION arctan# (y#, x#)
  37.  
  38. DEF fnyear% = VAL(RIGHT$(d$, 2))
  39. DEF fnhour% = VAL(LEFT$(t$, 2))
  40. DEF fnmonth% = VAL(LEFT$(d$, 2))
  41. DEF fnday% = VAL(MID$(d$, 4, 2))
  42. DEF fnminute% = VAL(MID$(t$, 4, 2))
  43. DEF fnsecond% = VAL(RIGHT$(t$, 2))
  44. DEF fnt$ (d%) = MID$(STR$(d% + 100), 3)
  45. DEF fnleap% (year%) = ((year% MOD 4) = 0)
  46. DEF fnradians (degrees) = degrees / 180# * pi#
  47. DEF fndegrees (radians) = radians / pi# * 180#
  48. DEF fnarcsin (sine) = ATN(sine / SQR(1 - sine * sine))
  49. ' Greenwich Mean Sidereal Time
  50. DEF fngmst# (YR%) = (99.6367# - .2387# * (YR% - 1989) + .9856# * INT((YR% - 1989) / 4#)) / 360!
  51.  
  52. CLEAR
  53.  
  54. maxsats% = 17
  55. DIM name$(maxsats%)
  56. DIM epochyear%(maxsats%), TLERevolution&(maxsats%), satl(maxsats%, 2)
  57. DIM epochjulianday#(maxsats%), inclination#(maxsats%), TLERAAN#(maxsats%), eccentricity#(maxsats%), TLEArgOfPerigee#(maxsats%)
  58. DIM TLEMA#(maxsats%), MeanMotion#(maxsats%)
  59. DIM BeaconFrequency&(maxsats%), a0(maxsats%)
  60. DIM sat%(7), pkt%(7), kep%(7) ' Image arrays (9 x 5 pixels)
  61.  
  62. DIM cc(3, 2) ' Used in orbit determination routines
  63.    
  64. ' ****** NUMERIC CONSTANTS ******
  65. pi# = 3.141592653589793# ' Value of PI
  66. r0# = 6378.16#           ' Earth's radius
  67. f# = 1# / 298.16#        ' 1/Earth flattening coef.
  68. g0# = 75369793000000#    ' GM of Earth in (orbits/day)^2/km^3
  69. g1# = 1.0027379093#      ' Sidereal/Solar time rate ratio
  70. c# = 299792.5#           ' Speed of light (km/sec).
  71.  
  72. ' String constants
  73. tlefile$ = "USAT.TLE"   ' Default NASA/NORAD element file
  74.  
  75. ' YOU WILL NEED THE FILES USATCGA.MAP, USAT.TLE AND USATGRND.DAT
  76. ' TO RUN THIS PROGRAM
  77. ' **** SATAUS Menuprogramm de OE1HSI VERSION 1.5  26.JAN. 1985
  78.  
  79. ' ****** ESTABLISH SATELLITE ELEMENT MATRIX ******
  80. PRINT "Reading satellite elements from "; tlefile$
  81. numsats% = readtle%(tlefile$, name$(), epochyear%(), epochjulianday#(), inclination#(), TLERAAN#(), eccentricity#(), TLEArgOfPerigee#(), TLEMA#(), MeanMotion#(), TLERevolution&())
  82.  
  83. ' Initialize beacon matrix and check age of satellite elements.
  84. FOR i% = 1 TO numsats%
  85.   BeaconFrequency&(i%) = 0
  86.   IF epochyear%(i%) < fnyear% - 1 THEN
  87.     PRINT "Warning:  Elements for satellite "; trim(name$(i%)); " are very old."
  88.     PRINT "You should update your "; tlefile$; " file."
  89.   END IF
  90. NEXT
  91.  
  92. ' Initialize the ground station data
  93. CALL readgrnd
  94.  
  95. ' ****** DRAW AND STORE SATELLITE SYMBOLS ******
  96. CLS
  97. SCREEN 2
  98. LINE (4, 0)-(4, 4)
  99. LINE (0, 2)-(8, 2)
  100. GET (0, 0)-(8, 4), sat%
  101. PUT (0, 0), sat%
  102. CLS
  103. LINE (4, 1)-(4, 3)
  104. LINE (3, 2)-(5, 2)
  105. GET (0, 0)-(8, 4), pkt%
  106. PUT (0, 0), pkt%
  107.  
  108. ' ***** Initializations done. *****
  109.  
  110. 50
  111. ' ****** DERIVED CONSTANTS ******
  112. ' These values depend on the ground station data.
  113. ' C$=GROUND STATION CALL+LOCATION STRING
  114. c$ = trim$(gsid$) + " " + trim$(gsloc$)
  115. pi2# = 2 * pi#
  116. SinGSLat# = SIN(fnradians(gslat!))    ' sin/cos of ground station latitude
  117. CosGSLat# = COS(fnradians(gslat!))
  118. SinGSLong# = SIN(fnradians(-gslong!))  ' sin/cos of ground station longitude
  119. CosGSLong# = COS(fnradians(gslong!))
  120. R9 = r0# * (1 - (f# / 2) + (f# / 2) * COS(2 * fnradians(gslat!))) + gsheight% / 1000
  121. L8 = ATN((1 - f#) ^ 2 * SinGSLat# / CosGSLat#)
  122. GSzcoord = R9 * SIN(L8)
  123. GSxcoord = R9 * COS(L8) * CosGSLong#
  124. GSycoord = R9 * COS(L8) * SinGSLong#
  125.  
  126. 80
  127. CLS
  128. SCREEN 0
  129. KEY(9) OFF
  130. KEY(10) OFF
  131. PRINT "                        USAT90.BAS - Version 1989-91"
  132. PRINT
  133. PRINT "======================================================================"
  134. PRINT
  135. PRINT "            SELECT ONE OF THE FOLLOWING OPTIONS:"
  136. PRINT
  137. PRINT " (P) ORBITAL PREDICTION PROGRAM"
  138. PRINT
  139. PRINT " (R) REALTIME TRACKING AND HIGH RESOLUTION SCREEN"
  140. PRINT
  141. PRINT " (S) Display ELEMENTS OF SATELLITES"
  142. PRINT
  143. PRINT " (G) CHANGE OR ENTER GROUNDSTATION DATA"
  144. PRINT
  145. PRINT " (D) RETURN TO DOS"
  146. PRINT
  147. PRINT
  148. PRINT "ENTER SELECTION (P,R,C,G,D)--> "
  149. CALL getkey(k$)
  150. ON INSTR("PRSGD", k$) GOTO 4820, 2090, 240, 1620, 9000
  151. BEEP
  152. GOTO 50
  153.  
  154. 240
  155. ' ****** SATFILE.BAS - VERSION 1.0, ISSUE 1.0 - HSIMODIF.1/25/85 ******
  156. CLS
  157. DO
  158.   PRINT "Elements of the following satellites are available:"
  159.   PRINT
  160.   FOR j% = 1 TO numsats%
  161.     PRINT name$(j%),
  162.   NEXT
  163.   PRINT
  164.   ' ****** FIND SATELLITE RECORD ******
  165.   DO
  166.     PRINT : INPUT "Which satellite"; v1$: IF v1$ = "" THEN 80
  167.     v1$ = UCASE$(trim$(v1$))
  168.     FOR j% = 1 TO numsats%
  169.       IF UCASE$(trim$(name$(j%))) = v1$ THEN
  170.         CLS
  171.         ' ****** DISPLAY SATELLITE RECORD ******
  172.         PRINT "Satellite       = "; name$(j%)
  173.         PRINT "Epoch year      = "; epochyear%(j%)
  174.         PRINT "Epoch day       = "; epochjulianday#(j%)
  175.         PRINT "Inclination     = "; inclination#(j%)
  176.         PRINT "R.A.A.N.        = "; TLERAAN#(j%)
  177.         PRINT "Eccentricity    = "; eccentricity#(j%)
  178.         PRINT "Arg. of perigee = "; TLEArgOfPerigee#(j%)
  179.         PRINT "Mean anomaly    = "; TLEMA#(j%)
  180.         PRINT "Mean motion     = "; MeanMotion#(j%)
  181.         PRINT "Epoch orbit no. = "; TLERevolution&(j%)
  182.         PRINT "Beacon freq.    = "
  183.         PRINT
  184.         EXIT DO
  185.       END IF
  186.     NEXT
  187.     PRINT "Satellite not found."
  188.   LOOP
  189. LOOP
  190.  
  191. 1620
  192. ' ******* Groundsation data change v.1.0 OE1HSI  jan.-1985**********
  193. CLS
  194. PRINT "CURRENT GROUND STATION DATA"
  195. PRINT
  196. GOSUB 2080
  197. DO
  198.   PRINT "Do you want to CHANGE this DATA ? (Y/N)"
  199.   CALL getkey(k$)
  200.   ON INSTR("YN", k$) GOTO 1780, 2060
  201.   BEEP
  202. LOOP
  203. 1710
  204. PRINT
  205. PRINT
  206. GOSUB 2080
  207. DO
  208.   PRINT "Do you want to make further changes ? (Y/N) "
  209.   CALL getkey(k$)
  210.   ON INSTR("YN", k$) GOTO 1780, 2050
  211.   BEEP
  212. LOOP
  213. 1780
  214. PRINT
  215. PRINT "ENTER NEW DATA OR <RETURN> FOR UNCHANGED DATA":
  216. PRINT
  217. PRINT "Ground Station Identifier                    : ";
  218. LINE INPUT u$
  219. IF trim$(u$) <> "" THEN gsid$ = trim$(u$)
  220. PRINT "Location of station                          : ";
  221. LINE INPUT u$
  222. IF trim$(u$) <> "" THEN gsloc$ = trim$(u$)
  223. PRINT
  224. PRINT "Longitude WEST of Greenwich (max +360) or East of Greenw. entered as -0 to -180"
  225. PRINT
  226. INPUT "Enter (with decimals)                 : ", u$
  227. IF trim$(u$) <> "" THEN gslong! = VAL(u$)
  228. IF gslong! < 0 THEN gslong! = 360 + gslong!
  229. PRINT
  230. PRINT "LATITUDE NORTH of Equator + (max 90) SOUTH of Equator - (max 90)"
  231. PRINT
  232. INPUT "ENTER (With decimals)                 : ", u$
  233. IF trim$(u$) <> "" THEN gslat! = VAL(u$)
  234. PRINT
  235. INPUT "Station height above sea level in meters     : ", u$
  236. IF trim$(u$) <> "" THEN gsheight% = VAL(u$)
  237. PRINT
  238. PRINT "Correction from local computer time to UTC"
  239. PRINT
  240. INPUT "Enter (in hours with decimals)        : ", u$
  241. IF trim$(u$) <> "" THEN tzcor! = VAL(u$)
  242. GOTO 1710
  243. 2050
  244. PRINT
  245. PRINT "DATA SAVED AS USATGRND.DAT"
  246. OPEN "USATGRND.DAT" FOR OUTPUT AS #1
  247. PRINT #1, gsid$
  248. PRINT #1, gsloc$
  249. PRINT #1, gslong!, gslat!, gsheight%, tzcor!
  250. CLOSE #1
  251. GOTO 50
  252. 2060
  253. PRINT
  254. PRINT "DATA NOT CHANGED"
  255. GOTO 50
  256.  
  257. 2080
  258.   PRINT "Ground Station Identifier            : "; gsid$
  259.   PRINT "Location name                        : "; gsloc$
  260.   PRINT USING "West longitude (degrees)             : +###.##"; gslong!
  261.   PRINT USING "Latitude (degrees)                   :  +##.##"; gslat!
  262.   PRINT USING "Height above Mean Sea Level (meters) :   #####"; gsheight%
  263.   PRINT USING "Local time correction to UTC (hours) :  +##.##"; tzcor!
  264.   PRINT
  265. RETURN
  266. '**** END PROGRAM GROUNDSTATION DATA CHANGE/STORAGE OE1HSI  JAN. 1985 ****
  267.  
  268. 2090
  269. ' ****** ORBITS2 - VERSION 1.0, ISSUE 1.2 -11/1/83 ******
  270. KEY OFF
  271. SCREEN 2, 0
  272. CLS
  273.  
  274. ' Initialize D$ and T$ from date$ and time$
  275. CALL timezone(DATE$, TIME$, d$, t$, tzcor!)
  276. ' ****** SET UP KEY TRAPPING ******
  277. ON KEY(9) GOSUB 4760
  278. KEY(9) STOP
  279. ON KEY(10) GOSUB 4790
  280. KEY(10) OFF
  281. flg9 = 0
  282. flg10 = 0
  283. GOSUB 4290
  284.    
  285. ' ****** ORBIT DETERMINATION LOOP STARTS HERE ******
  286. q$ = ""
  287. DO UNTIL q$ = CHR$(27)
  288.   FOR j% = 1 TO numsats%
  289.     q$ = UCASE$(INKEY$)
  290.     IF q$ = CHR$(27) THEN EXIT FOR
  291.     CALL timezone(DATE$, TIME$, d$, t$, tzcor!)
  292.     t = INT(30.55 * (fnmonth% + 2)) - 2 * (INT(.1 * (fnmonth% + 7))) - 91
  293.     IF fnmonth% > 2 AND fnleap%(fnyear%) THEN t = t + 1
  294.     IF epochyear%(j%) = fnyear% - 1 THEN
  295.       t = t + 365
  296.       IF fnleap%(epochyear%(j%)) THEN t = t + 1
  297.     END IF
  298.     t = t + fnday% + fnhour% / 24# + fnminute% / 1440# + fnsecond% / 86400#
  299.     GOSUB 6780
  300.     GOSUB 6910
  301.     IF flg9 THEN CALL plotloc(longitude!, latitude!, plotx%, ploty%): GOSUB 3890
  302.     GOSUB 4040
  303.   NEXT
  304. LOOP
  305. GOTO 50
  306.    
  307. ' ****** PUT SATELLITE ON SCREEN ROUTINE ******
  308. 3890
  309.   GET (plotx% - 4, ploty% - 2)-(plotx% + 4, ploty% + 2), kep%
  310.   PUT (plotx% - 4, ploty% - 2), sat%, PRESET
  311.   PUT (plotx% - 4, ploty% - 2), sat%
  312.   PUT (plotx% - 4, ploty% - 2), kep%, PSET
  313.   PUT (plotx% - 4, ploty% - 2), sat%
  314.   IF FLG0 <> 0 THEN
  315.     PUT (satl(j%, 1), satl(j%, 2)), sat%
  316.     PUT (satl(j%, 1), satl(j%, 2)), pkt%, OR
  317.   END IF
  318.   satl(j%, 1) = plotx% - 4
  319.   satl(j%, 2) = ploty% - 2
  320.   IF j% = numsats% THEN FLG0 = 1
  321. RETURN
  322.  
  323. ' ****** PRINT SATELLITE DETAILS ROUTINE ******
  324. 4040
  325.   KEY(9) ON
  326.   KEY(10) ON
  327. 4050
  328.   KEY(10) STOP
  329.   KEY(9) STOP
  330.   IF FLGK = 1 THEN GOSUB 4270
  331.   IF flg9 = 0 THEN
  332.     LOCATE 2, 49
  333.     PRINT d$
  334.     LOCATE 2, 63
  335.     PRINT t$
  336.     IF elevation! >= 0 THEN
  337.       IF elevation! > 0 AND elevation! < 1 THEN BEEP
  338.     END IF
  339.     IF j% + 6 > 23 GOTO 4230
  340.     LOCATE j% + 6, 15
  341.     PRINT SPACE$(50)
  342.     LOCATE j% + 6, 8 - (LEN(trim$(name$(j%))) / 2)
  343.     PRINT trim$(name$(j%))
  344.     LOCATE j% + 6, 15
  345.   ELSE
  346.     IF flg10 = 1 THEN GOSUB 4540
  347.     IF flg10 <> 0 THEN
  348.       IF UCASE$(trim$(name$(j%))) <> u$ THEN RETURN
  349.       LOCATE 25, 69
  350.       PRINT "SELECTED";
  351.     END IF
  352.     LOCATE 25, 1
  353.     PRINT SPACE$(68);
  354.     LOCATE 25, 8 - (LEN(trim$(name$(j%))) / 2)
  355.     PRINT trim$(name$(j%));
  356.     LOCATE 25, 15
  357.   END IF
  358.   PRINT USING "###   ###  #####    #####  ###.#   ###.#    ######"; azimuth!; elevation!; range#; (r - r0#); latitude!; longitude!; revolution&;
  359. 4230
  360.   IF flg9 <> 0 THEN LOCATE 20, 48: PRINT t$;
  361. RETURN
  362.    
  363. ' ****** SET UP SCREEN DISPLAY ROUTINE ******
  364. 4270
  365.   CLS
  366.   FLG0 = 0
  367.   FLGK = 0
  368. 4290
  369.   IF flg9 = 1 THEN
  370.     SCREEN 2, 0
  371.     DEF SEG = &HB800
  372.     BLOAD "USATCGA.MAP", 0
  373.     DEF SEG
  374.     CALL plotloc(gslong!, gslat!, plotx%, ploty%)
  375.     CIRCLE (plotx%, ploty%), 2
  376.     LOCATE , , 1, 12, 13
  377.     L1 = 22
  378.     L2 = 23
  379.     L3 = 24
  380.     LOCATE 20, 3
  381.     PRINT "Data for Groundstation             At Time=           UTC  On: "; d$
  382.     LOCATE 20, 26
  383.     PRINT gsid$
  384.   ELSE
  385.     ' FLG9=0
  386.     SCREEN 0, 1
  387.     LOCATE 1, 40 - LEN(c$) / 2, 0
  388.     PRINT c$
  389.     LOCATE 2, 5
  390.     PRINT "Real-time satellite tracking coordinates on            at          UTC"
  391.     LOCATE 25, 3
  392.     PRINT "F9 TOGGLES THE GRAPH/TABLE     F10 TO SELECT SINGLE SAT IN GRAPH   ESC TO END";
  393.     L1 = 4
  394.     L2 = 5
  395.     L3 = 6
  396.   END IF
  397.   LOCATE L1, 3
  398.   PRINT "  NAME OR   AZ    EL   RANGE   HEIGHT   LAT.   LONG.    ORBIT"
  399.   LOCATE L2, 3
  400.   PRINT "DESIGNATOR  DEG   DEG    KM      KM     DEG     DEG       NO."
  401.   LOCATE L3, 3
  402.   PRINT "----------  ---   ---  -----   ------  -----   -----    ------";
  403. RETURN
  404.    
  405. ' ****** SELECT SINGLE SATELLITE ROUTINE ******
  406. 4540
  407.   LOCATE 25, 1
  408.   PRINT SPACE$(79);
  409.   LOCATE 25, 1
  410.   INPUT ; "WHICH SATELLITE? (CR FOR ALL)"; i$
  411.   i$ = UCASE$(trim$(i$))
  412.   IF i$ = "" THEN
  413.     flg10 = 0
  414.   ELSE
  415.     u$ = i$
  416.     flg10 = 2
  417.   END IF
  418.   LOCATE 25, 1
  419.   PRINT SPACE$(79);
  420. RETURN
  421.  
  422. 4760 ' toggle FLG9
  423.   flg9 = -(flg9 = 0)
  424.   FLGK = 1
  425. RETURN 4050
  426.  
  427. 4790 flg10 = flg9
  428. RETURN 4050
  429.  
  430. 4820
  431.   '****** ORBIT2 - VERSION 2.0, ISSUE 1.0/HSI - 17/01/85 *****
  432.   KEY OFF
  433.   SCREEN 0
  434.   CLS
  435.   ' ****** HOUSEKEEPING ITEMS ******
  436.   C8$ = CHR$(10) + CHR$(10) + CHR$(10) + CHR$(10)
  437.   C9$ = CHR$(12) + CHR$(7)
  438.  
  439.   DO
  440.     PRINT "Enter start date (mm-dd-yy):  ";
  441.     u$ = gets$
  442.     year% = VAL(RIGHT$(u$, 2))
  443.     month% = VAL(LEFT$(u$, 2))
  444.     day% = VAL(MID$(u$, 4, 2))
  445.     IF month% >= 1 AND month% <= 12 AND day% >= 1 AND day% <= 31 AND year% >= 91 AND MID$(u$, 3, 1) = "-" AND MID$(u$, 6, 1) = "-" THEN EXIT DO
  446.   LOOP
  447.   YY = year% + 1900
  448.   t$ = fnt$(year%) + "/" + fnt$(month%) + "/" + fnt$(day%) + " at "
  449.   D8 = day% + INT(30.55 * (month% + 2)) - 2 * (INT(.1 * (month% + 7))) - 91
  450.   IF month% > 2 THEN IF fnleap%(year%) THEN D8 = D8 + 1
  451.   PRINT "    Day #"; D8
  452.   PRINT
  453.   PRINT "Enter start time (hh:mm):  ";
  454.   DO
  455.     u$ = gets$
  456.     h% = VAL(LEFT$(u$, 2))
  457.     m% = VAL(RIGHT$(u$, 2))
  458.     IF h% >= 0 AND h% <= 23 AND m% >= 0 AND m% <= 59 AND MID$(u$, 3, 1) = ":" THEN EXIT DO
  459.     PRINT "Valid format is:  hh:mm"
  460.     PRINT "Re-enter:  ";
  461.   LOOP
  462.   T7 = D8 + h% / 24# + m% / 1440#
  463.   t$ = t$ + fnt$(h%) + fnt$(m%) + ":00 H"
  464.   PRINT "Enter duration in hours and minutes (hh:mm):  ";
  465.   DO
  466.     u$ = gets$
  467.     h% = VAL(LEFT$(u$, 2))
  468.     m% = VAL(RIGHT$(u$, 2))
  469.     IF h% >= 0 AND m% >= 0 AND m% <= 59 AND MID$(u$, 3, 1) = ":" THEN EXIT DO
  470.     PRINT "Valid format is:  hh:mm"
  471.     PRINT "Re-enter:  ";
  472.   LOOP
  473.   T8 = T7 + h% / 24# + m% / 1440#
  474.   PRINT USING "    From ###.####### to ###.#######"; T7; T8
  475.   DO
  476.     PRINT "Enter time increment (minutes):  ";
  477.     INPUT m%
  478.     IF m% > 0 THEN EXIT DO
  479.   LOOP
  480.   T9 = m% / 1440#
  481.   PRINT
  482.   INPUT "MINIMUM ELEVATION ? (DEFAULT 0) Deg. = ", E8
  483.   DO
  484.     PRINT
  485.     INPUT "Output to Printer (P) or Screen (S) ?-->", P$
  486.     P$ = UCASE$(P$)
  487.     IF (P$ = "P" OR P$ = "S") THEN
  488.       EXIT DO
  489.     ELSE BEEP
  490.     END IF
  491.   LOOP
  492.   IF P$ = "P" THEN outfile$ = "prn:" ELSE outfile$ = "scrn:"
  493.   CLOSE #2
  494.   OPEN outfile$ FOR OUTPUT AS #2
  495.   IF outfile$ = "prn:" THEN
  496.     PRINT "Put the printer on-line and align the page,"
  497.     PRINT "then press any key when ready."
  498.     DO
  499.       getkey (u$)
  500.       IF u$ <> "" THEN EXIT DO
  501.     LOOP
  502.   END IF
  503.  
  504. ' ****** SELECT SATELLITE FROM MENU ******
  505.   CLS
  506.   PRINT "SATELLITE SELECTION MENU"
  507.   PRINT
  508.   FOR j% = 1 TO numsats%
  509.     PRINT name$(j%),
  510.   NEXT
  511.   PRINT
  512.   DO
  513.     PRINT
  514.     INPUT "SELECT satellite >", y3$
  515.     y3$ = UCASE$(trim$(y3$))
  516.     FOR j% = 1 TO numsats%
  517.       IF UCASE$(trim$(name$(j%))) = y3$ THEN EXIT DO
  518.     NEXT
  519.     PRINT "Satellite not found."
  520.   LOOP
  521.   PRINT
  522.   PRINT "Doppler calculated for frequ. = "; BeaconFrequency&(j%); " MHz"
  523.   INPUT " Change frequency to: (0 for default) ", d
  524.   IF d <> 0 THEN BeaconFrequency&(j%) = d
  525.   PRINT #2,
  526.   PRINT #2, "Orbital ELEMENTS for "; name$(j%)
  527.   PRINT #2,
  528.   PRINT #2, "Reference epoch = "; epochyear%(j%); " +"; epochjulianday#(j%)
  529.   PRINT #2, "Starting  epoch = "; year%; " +"; T7; " = "; t$
  530.   PRINT #2,
  531.   PRINT #2, "Parameter"; TAB(20); "Reference"; TAB(40); "Starting"
  532.   t = T7
  533.   IF epochyear%(j%) = year% - 1 THEN
  534.     t = t + 365
  535.     T8 = T8 + 365
  536.     IF fnleap%(epochyear%(j%)) THEN t = t + 1: T8 = T8 + 1
  537.   END IF
  538.   ' Calculate start time parameters
  539.   GOSUB 6780
  540.   PRINT #2, "Orbit Number "; TAB(20); TLERevolution&(j%); TAB(40); revolution&
  541.   PRINT #2, "Mean Anomaly "; TAB(20); TLEMA#(j%); TAB(40); fndegrees(MARadians#)
  542.   PRINT #2, "Inclination  "; TAB(20); inclination#(j%)
  543.   PRINT #2, "Eccentricity "; TAB(20); eccentricity#(j%)
  544.   PRINT #2, "Mean Motion  "; TAB(20); MeanMotion#(j%)
  545.   PRINT #2, "S.M.A.,km    "; TAB(20); a0(j%)
  546.   PRINT #2, "Arg. Perigee "; TAB(20); TLEArgOfPerigee#(j%); TAB(40); ArgOfPerigee#
  547.   PRINT #2, "R. A. A. N.  "; TAB(20); TLERAAN#(j%); TAB(40); RAAN#
  548.   PRINT #2, "Freq.,MHz    "; TAB(20); BeaconFrequency&(j%)
  549.   OldRevolution& = -1
  550.   OldJulianDay% = -1
  551.   
  552.   '****** COMPUTATION LOOP ******
  553. 6400
  554.   t = t + T9
  555.   JulianDay% = INT(t)
  556.   IF JulianDay% <> OldJulianDay% THEN
  557.     OldJulianDay% = -1
  558.     OldRevolution& = -1
  559.   END IF
  560.   GOSUB 6780
  561.   GOSUB 6910
  562.   IF elevation! >= E8 THEN
  563.     IF JulianDay% <> OldJulianDay% OR revolution& <> OldRevolution& THEN
  564.       IF JulianDay% <> OldJulianDay% THEN
  565.         GOSUB 6690
  566.         OldJulianDay% = JulianDay%
  567.         PRINT #2, " U.T.C.    AZ    EL  DOPPLER   RANGE   HEIGHT  LAT.  LONG.  PHASE"
  568.         PRINT #2, "HHMM:SS   deg   deg    Hz       km       km    deg    deg   <256>"
  569.       END IF
  570.       PRINT #2, TAB(21); "- - - ORBIT #"; revolution&; "- - -"
  571.     END IF
  572.     OldRevolution& = revolution&
  573.     T4 = t - JulianDay%
  574.     S4 = INT(T4 * 86400)
  575.     H4 = INT(S4 / 3600 + .000001)
  576.     M4 = INT((S4 - H4 * 3600) / 60 + .000001)
  577.     S4 = S4 - 3600 * H4 - 60 * M4
  578.     doppler# = -BeaconFrequency&(j%) * 1000000 * R8 / c#
  579.     PRINT #2, fnt$(H4) + fnt$(M4) + ":" + fnt$(S4);
  580.     PRINT #2, USING "   ###   ###  #####"; azimuth!; elevation!; doppler#;
  581.     PRINT #2, USING "    #####    #####  ###.#  ###.#  ###"; range#; (r - r0#); latitude!; longitude!; phase%
  582.   END IF
  583.   IF t < T8 THEN GOTO 6400
  584.   PRINT #2, C9$
  585.   PRINT "Do YOU have another INQUIRY  ?  (Y/N) "
  586.   PRINT
  587.   PRINT "Else you return to the MAIN MENU !"
  588.   PRINT
  589.   DO
  590.     CALL getkey(k$)
  591.     ON INSTR("YN", k$) GOTO 4820, 6670
  592.     BEEP
  593.   LOOP
  594. 6670
  595.   CLOSE #2
  596. GOTO 80
  597.     
  598. '****** PAGE HEADER SUBROUTINE ******
  599. 6690
  600.   PRINT #2, C9$; c$; "   Lat.="; gslat!; "  W.Long.="; gslong!; "  Ht.="; gsheight%;
  601.   pagenum% = pagenum% + 1
  602.   PRINT #2, TAB(70); "Page # "; pagenum%
  603.   PRINT #2, TAB(15); " - - - Minimum Elevation = "; E8; "Deg. - - -"
  604.   PRINT #2,
  605.   PRINT #2, TAB(14); "- - - DAY #"; JulianDay%; "- - -  "; juliantodate$(JulianDay%, INT(YY)); "  - - -"
  606.   PRINT #2,
  607. RETURN
  608.  
  609. '****** ORBIT DETERMINATION AND UTILITY ROUTINES ******
  610.  
  611. 6780
  612.   a0(j%) = ((g0# / (MeanMotion#(j%) * MeanMotion#(j%))) ^ (1 / 3))
  613.   e2 = 1 - eccentricity#(j%) * eccentricity#(j%)
  614.   E1 = SQR(e2)
  615.   TLEOrbitPos# = TLEMA#(j%) / 360 + TLERevolution&(j%)
  616.   K2 = 9.95 * ((r0# / a0(j%)) ^ 3.5) / (e2 * e2)
  617.   S1 = SIN(fnradians(inclination#(j%)))
  618.   C1 = COS(fnradians(inclination#(j%)))
  619.   RAAN# = TLERAAN#(j%) - (t - epochjulianday#(j%)) * K2 * C1
  620.   S0 = SIN(fnradians(RAAN#))
  621.   C0 = COS(fnradians(RAAN#))
  622.   ArgOfPerigee# = TLEArgOfPerigee#(j%) + (t - epochjulianday#(j%)) * K2 * (2.5 * (C1 * C1) - .5)
  623.   S2 = SIN(fnradians(ArgOfPerigee#))
  624.   C2 = COS(fnradians(ArgOfPerigee#))
  625.   cc(1, 1) = (C2 * C0) - (S2 * S0 * C1)
  626.   cc(1, 2) = -(S2 * C0) - (C2 * S0 * C1)
  627.   cc(2, 1) = (C2 * S0) + (S2 * C0 * C1)
  628.   cc(2, 2) = -(S2 * S0) + (C2 * C0 * C1)
  629.   cc(3, 1) = (S2 * S1)
  630.   cc(3, 2) = (C2 * S1)
  631.   OrbitPos# = TLEOrbitPos# + MeanMotion#(j%) * (t - epochjulianday#(j%))
  632.   revolution& = INT(OrbitPos#)
  633.   phase% = INT((OrbitPos# - revolution&) * 256)
  634.   MARadians# = (OrbitPos# - revolution&) * pi2#
  635. RETURN
  636.  
  637. 6910
  638.   E = MARadians# + eccentricity#(j%) * SIN(MARadians#) + .5 * (eccentricity#(j%) * eccentricity#(j%)) * SIN(2 * MARadians#)
  639.   DO
  640.     S3 = SIN(E)
  641.     C3 = COS(E)
  642.     R3 = 1 - eccentricity#(j%) * C3
  643.     M1 = E - eccentricity#(j%) * S3
  644.     M5 = M1 - MARadians#
  645.     IF ABS(M5) < .000001 THEN EXIT DO ELSE E = E - M5 / R3
  646.   LOOP
  647.   X0 = a0(j%) * (C3 - eccentricity#(j%))
  648.   y0 = a0(j%) * E1 * S3
  649.   r = a0(j%) * R3
  650.   X1 = X0 * cc(1, 1) + y0 * cc(1, 2)
  651.   y1 = X0 * cc(2, 1) + y0 * cc(2, 2)
  652.   Z1 = X0 * cc(3, 1) + y0 * cc(3, 2)
  653.   G7 = t * g1# + fngmst#(1900 + epochyear%(j%))
  654.   G7 = (G7 - INT(G7)) * pi2#
  655.   S7 = -SIN(G7)
  656.   C7 = COS(G7)
  657.   x = (X1 * C7) - (y1 * S7)
  658.   y = (X1 * S7) + (y1 * C7)
  659.   z = Z1
  660.   X5 = (x - GSxcoord)
  661.   Y5 = (y - GSycoord)
  662.   Z5 = (z - GSzcoord)
  663.   range# = SQR(X5 * X5 + Y5 * Y5 + Z5 * Z5)
  664.   ' Speed to/from ground station
  665.   IF OldT# <> t THEN R8 = ((range# - OldRange#) / (t - OldT#)) / 86400 ELSE R8 = -900000000
  666.   OldRange# = range#
  667.   OldT# = t
  668.   ' Translate to ground station coordinate system
  669.   z8 = (X5 * CosGSLong# * CosGSLat#) + (Y5 * SinGSLong# * CosGSLat#) + (Z5 * SinGSLat#)
  670.   x8 = -(X5 * CosGSLong# * SinGSLat#) - (Y5 * SinGSLong# * SinGSLat#) + (Z5 * CosGSLat#)
  671.   y8 = (Y5 * CosGSLong#) - (X5 * SinGSLong#)
  672.   elevation! = fndegrees(fnarcsin(z8 / range#))
  673.   azimuth! = fndegrees(arctan(y8, x8))
  674.   longitude! = 360 - fndegrees(arctan(y, x))
  675.   latitude! = fndegrees(fnarcsin(z / r))
  676. RETURN
  677.  
  678. 9000 END
  679.  
  680. FUNCTION arctan (y, x)
  681. ' Modified arctangent function.  Returns a value between 0 and 2 PI.
  682. SHARED pi#, pi2#
  683.  
  684. IF x > 0 AND y >= 0 THEN
  685.   arctan = ATN(y / x)
  686. ELSEIF x > 0 AND y < 0 THEN
  687.   arctan = ATN(y / x) + pi2#
  688. ELSEIF x < 0 THEN
  689.   arctan = ATN(y / x) + pi#
  690. ELSEIF y >= 0 THEN
  691.   arctan = pi# / 2
  692. ELSE arctan = 1.5# * pi#
  693. END IF
  694. END FUNCTION
  695.  
  696. DEFSNG A-Z
  697. ' Converts a date string to the julian day of the year
  698. FUNCTION datetojulian% (d$)
  699.  
  700. day% = VAL(MID$(d$, 4, 2))
  701.  
  702. SELECT CASE LEFT$(d$, 2)
  703. CASE "01"
  704. CASE "02"
  705.   day% = day% + 31
  706. CASE "03"
  707.   days% = day% + 59
  708. CASE "04"
  709.   day% = day% + 90
  710. CASE "05"
  711.   day% = day% + 120
  712. CASE "06"
  713.   day% = day% + 151
  714. CASE "07"
  715.   day% = day% + 181
  716. CASE "08"
  717.   day% = day% + 212
  718. CASE "09"
  719.   day% = day% + 243
  720. CASE "10"
  721.   day% = day% + 273
  722. CASE "11"
  723.   day% = day% + 304
  724. CASE "12"
  725.   day% = day% + 334
  726. END SELECT
  727.  
  728. IF fnleap%(VAL(RIGHT$(d$, 4))) AND LEFT$(d$, 2) > "02" THEN
  729.   datetojulian% = day% + 1
  730. ELSE
  731.   datetojulian% = day%
  732. END IF
  733.  
  734. END FUNCTION
  735.  
  736. SUB getkey (k$)
  737.  
  738. k$ = ""
  739. WHILE k$ = ""
  740.         k$ = UCASE$(INKEY$)
  741. WEND
  742.  
  743. END SUB
  744.  
  745. ' Gets an input line, trims leading and trailing blanks
  746. ' then converts to upper-case before returning.
  747. FUNCTION gets$
  748. LINE INPUT i$
  749. gets$ = UCASE$(trim$(i$))
  750. END FUNCTION
  751.  
  752. ' Converts a julian day and year to a date string.
  753. FUNCTION juliantodate$ (day%, year%)
  754.  
  755. IF day% = 60 THEN
  756.   IF fnleap%(year%) THEN tempdate$ = "02-29" ELSE tempdate$ = "03-01"
  757.   GOTO appendyear
  758. ELSE IF fnleap%(year%) AND day% > 60 THEN day% = day% - 1
  759. END IF
  760.  
  761. SELECT CASE day%
  762. CASE 1 TO 31
  763.   tempdate$ = "01-" + MID$(STR$(day% + 100), 3)
  764. CASE 32 TO 59
  765.   tempdate$ = "02-" + MID$(STR$(day% - 31 + 100), 3)
  766. CASE 60 TO 90
  767.   tempdate$ = "03-" + MID$(STR$(day% - 59 + 100), 3)
  768. CASE 91 TO 120
  769.   tempdate$ = "04-" + MID$(STR$(day% - 90 + 100), 3)
  770. CASE 121 TO 151
  771.   tempdate$ = "05-" + MID$(STR$(day% - 120 + 100), 3)
  772. CASE 152 TO 181
  773.   tempdate$ = "06-" + MID$(STR$(day% - 151 + 100), 3)
  774. CASE 182 TO 212
  775.   tempdate$ = "07-" + MID$(STR$(day% - 181 + 100), 3)
  776. CASE 213 TO 243
  777.   tempdate$ = "08-" + MID$(STR$(day% - 212 + 100), 3)
  778. CASE 244 TO 273
  779.   tempdate$ = "09-" + MID$(STR$(day% - 243 + 100), 3)
  780. CASE 274 TO 304
  781.   tempdate$ = "10-" + MID$(STR$(day% - 273 + 100), 3)
  782. CASE 305 TO 334
  783.   tempdate$ = "11-" + MID$(STR$(day% - 304 + 100), 3)
  784. CASE 335 TO 365
  785.   tempdate$ = "12-" + MID$(STR$(day% - 334 + 100), 3)
  786. END SELECT
  787.  
  788. appendyear: juliantodate$ = tempdate$ + "-" + MID$(STR$(year%), 2)
  789.  
  790. END FUNCTION
  791.  
  792. SUB plotloc (longitude!, latitude!, xcoord%, ycoord%)
  793. ' Returns the screen x and y coordinates for a given map location
  794.  
  795. ycoord% = CINT(.7111 * (90 - latitude!) + 3)
  796. IF longitude! >= 0 AND longitude! <= 270 THEN
  797.   xcoord% = CINT(477 - longitude! * 1.7444)
  798. ELSE xcoord% = CINT(1105 - longitude! * 1.7444)
  799. END IF
  800.  
  801. END SUB
  802.  
  803. SUB readgrnd
  804. SHARED gsid$, gsloc$, gslong!, gslat!, gsheight%, tzcor!
  805.  
  806. ' Format of USATGRND.DAT:
  807. '   gsid$    = Station Identifier
  808. '   gsloc$   = Station location (name)
  809. '   gslong!   = Longitude in degrees
  810. '   gslat!    = Latitude in degrees
  811. '   gsheight% = Height above sea level in meters
  812. '   tzcor!    = Time zone correction to UTC
  813. OPEN "USATGRND.DAT" FOR INPUT AS #1
  814. LINE INPUT #1, gsid$
  815. LINE INPUT #1, gsloc$
  816. INPUT #1, gslong!, gslat!, gsheight%, tzcor!
  817. CLOSE #1
  818.  
  819. END SUB
  820.  
  821. FUNCTION readtle% (tlefile$, name$(), epochyear%(), epochjulianday#(), inclination#(), TLERAAN#(), eccentricity#(), TLEArgOfPerigee#(), TLEMA#(), MeanMotion#(), TLERevolution&())
  822. ' ****** ESTABLISH SATELLITE ELEMENT MATRIX ******
  823. ' tlefile$ conforms to the NASA/NORAD three-line element set format.
  824.  
  825. OPEN tlefile$ FOR INPUT AS #1
  826.  
  827. FOR i% = 1 TO UBOUND(name$)
  828.   IF EOF(1) THEN EXIT FOR
  829.   LINE INPUT #1, name$(i%)
  830.   IF trim$(name$(i%)) = "" THEN EXIT FOR
  831.   LINE INPUT #1, y3$
  832.   LINE INPUT #1, t0$
  833.   epochyear%(i%) = VAL(MID$(y3$, 19, 2))
  834.   epochjulianday#(i%) = VAL(MID$(y3$, 21, 12))
  835.   inclination#(i%) = VAL(MID$(t0$, 9, 8))
  836.   TLERAAN#(i%) = VAL(MID$(t0$, 18, 8))
  837.   eccentricity#(i%) = VAL("." + MID$(t0$, 27, 7))
  838.   TLEArgOfPerigee#(i%) = VAL(MID$(t0$, 35, 8))
  839.   TLEMA#(i%) = VAL(MID$(t0$, 44, 8))
  840.   MeanMotion#(i%) = VAL(MID$(t0$, 53, 11))
  841.   TLERevolution&(i%) = VAL(MID$(t0$, 64, 5))
  842. NEXT
  843.  
  844. CLOSE #1
  845. readtle% = i% - 1
  846.            
  847. END FUNCTION
  848.  
  849. ' Converts an input date and time string from one time zone to another.
  850. SUB timezone (d0$, t0$, d1$, t1$, tzcor!)
  851.  
  852. JulianDay% = datetojulian(d0$)
  853.  
  854. year% = VAL(RIGHT$(d0$, 4))
  855.  
  856. hour! = VAL(LEFT$(t0$, 2)) + INT(tzcor!) + (VAL(MID$(t0$, 4, 2)) + tzcor! - INT(tzcor!)) / 60
  857. SELECT CASE hour!
  858. CASE IS < 0
  859.   hour! = hour! + 24
  860.   JulianDay% = JulianDay% - 1
  861. CASE IS >= 24
  862.   hour! = hour! - 24
  863.   JulianDay% = JulianDay% + 1
  864. END SELECT
  865. t1$ = fnt$(INT(hour!)) + ":" + fnt$(INT(60 * (hour! - INT(hour!)))) + MID$(t0$, 6)
  866.  
  867. SELECT CASE JulianDay%
  868. CASE 0
  869.   year% = year% - 1
  870.   IF fnleap%(year%) THEN JulianDay% = 366 ELSE JulianDay% = 365
  871. CASE IS > 365
  872.   IF NOT fnleap%(year%) THEN JulianDay% = 1: year% = year% + 1
  873. END SELECT
  874.  
  875. d1$ = juliantodate$(JulianDay%, year%)
  876. END SUB
  877.  
  878. ' Trims leading and trailing spaces from a string.
  879. FUNCTION trim$ (i$)
  880.           
  881. trim$ = RTRIM$(LTRIM$(i$))
  882.  
  883. END FUNCTION
  884.  
  885.